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Q I Abstract 

We study the kinetic mean-field limits of the discrete systems of interacting particles 
used for halftoning of images in the sense of continuous-domain quantization. Under mild 
assumptions on the regularity of the interacting kernels we provide a rigorous derivation 
of the mean-field kinetic equation. Moreover, we study the energy of the system, show 
that it is a Lyapunov functional and prove that in the long time limit the solution tends 



in 

< 



, to an equilibrium given by a local minimum of the energy. In a special case we prove 



X 



that the equilibrium is unique and is identical to the prescribed image profile. This proves 
the consistency of the particle halftoning method when the number of particles tends to 
infinity. 
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■ 1 Introduction 

A halftoning method places black dots in an image in such a way that their density gives the 
impression of tone. For an illustration see Fig. [TJ Due to its various applications, halftoning 
is an active field of research and we refer to the recent papers [12^ |4] for deterministic and, 
^ ■ resp., stochastic point distributions. 

^ ■ In this paper, we consider a continuous-domain quantization method based on electrostatic- 

like principles studied in |14j . where the basic idea goes back to [12]. The method consists in 
considering a system of N particles with electrostatic-like interaction (repulsion) and exposed 
to an attractive external potential w with compact support in R*^ which represents the image 
to be approximated by points. In [l3| the authors considered the discrete energy functional 



N ^ ^ N N 



E{p)-=^l w{x)\pk- x\dx -^S2^\Pk-Pil (1) 



k=l k=l 1=1 
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Figure 1: Left: Original 256 x 256 image. Right: Halftoning result with m = 30150 points 
using the technique from |14] . 

where p ■= {pi, ■ ■ ■ ,Pn)^ ^ M.'^^'^, the Euclidean distance is denoted by | • |, 



and defined the evolution of the particle system as the corresponding gradient flow dtp € —dpE. 
They showed that for d = \ the energy functional is continuous and coercive, and calculated 
explicitly its minimizers. In the two-dimensional setting it is not possible to obtain explicit 
expressions for the minimizers anymore. Instead, the authors employed a difference of convex 
functions (DC) algorithm together with fast summation methods for non-equispaced knots to 
obtain a local minimum of the variational problem numerically. In [8] the function w was also 
considered on other sets such as T'^ or and kernels other than the Euclidian distance were 
used. This generalized approach stems from the study of optimal quadrature error functionals 
on reproducing kernel Hilbert spaces with respect to the quadrature knots and is closely related 
to so-called discrepancy functionals [7]. 

In this paper we study the kinetic mean-field limits of discrete systems like ([1]), that are ob- 
tained as the number of particles tends to infinity. First, in Section[2l we specify the interaction 
kernels which we will consider and provide some particular examples. Moreover, we introduce 
a generalized setting with different kernels for the attractive and repulsive interactions, and 
show that this new setting can be reduced to the previous one with appropriately modified 
data. In Section [3] we provide, under mild regularity conditions on the interaction kernels, a 
rigorous derivation of the mean field kinetic equation obtained in the limit as the number of 
particles tends to infinity. We show that the corresponding energy functional, obtained as a 
formal limit of the discrete energy, is a Lyapunov functional for the kinetic equation and that 
in the long-time limit the solution tends to an equilibrium, which is a local minimum of the 
energy. For a special choice of the interaction kernel, we are able to show that the equilibrium 
coincides with the prescribed image profile. This proves the consistency of the discrete halfton- 
ing method when the number of particles tends to infinity. Finally, in Section |4l we provide 
numerical examples confirming our theoretical results, and showing the behavior of the model 
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subject to different competitive attraction and repulsion terms as well as the consistency of 
the particle system with respect to its continuous kinetic limit. 



2 Energy Functionals with Kernels 

In the following, let Q he a domain, which we will be either the Euclidean space Q = 
or the torus Q = T'^, d > 1. A symmetric function iiT : Q x Q — )■ M is said to be positive 
semi-definite if for any € N points xi, . . . , x^r G Q and any a G \ {0} the relationship 
{K{xi,Xj))fj^ia > holds true, and positive definite if we have strict inequality. Let K : 
O X 17 — > M be a symmetric, positive semi-definite function, and Hk be the reproducing kernel 
Hilbert space (RKHS) associated with K, see jl]. Then we are interested in the functional 

N ^ N 

Ek{p) ■■= V / wix)Kip„x) dx - - V K{pi,pj) . (3) 



In the case = M we suppose that the function compactly supported. In [8] 

it was shown that this functional is related to the optimality of a certain quadrature rule for 
functions in Hk depending on the knots pi, i = 1, . . . , N . By the following remark, which was 
proved in [8j , slight modifications of the kernel do not change the minimizers of the functional 
Ek- 

Remark 1 Let K : Qxi} be a symmetric function and K{x, y) := aK{x, y) + b{K{x, 0) + 
K{0,y)) + c with a > and 6, c € M. Then the minimizers of Ek and Ej^ coincide. 

In this paper, we consider radial kernels on = M.'^, i.e., 

K{x,y) = ^{x-y) = <l>i\x-y\) (4) 

with <I> : [0, oo) — )• M. For the 1-periodic setting = T"^ we use the same notation, where 
|x — y| has to be replaced by the "periodic" distance min{|x — y|,l — |x — y|}. In the case 
= T'^ we consider tensor products of radial kernels. We call ip positive semi-definite (resp. 
positive definite) if the corresponding kernel is positive semi-definite (resp. positive definite). 
The functional of interest becomes then 

N ^ N 

EAp) '■=y2 ^i^MPi -x) dx--y] ip{pi - Pj) . (5) 

i=l ^ i,j=l 

Example. We give some interesting examples of positive semi-definite kernels, see [TB] and 

m. 



1. Let O := M . Then the functions 

V^(x) = (l-lx|);, r> 



+ 1, 



fix) = (e^ + |xp) ^ , ^ ^ 2 (inverse multiquadrics) 
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are positive definite. Next, consider the conditionally positive definite radial kernels of 
order 1 defined by 

ip{x) := -Ixl"", < r < 2, 

:= —{e^ + \x^y , < T < 1, (multiquadrics) . 

The kernels K{x, y) = (p{x — y) are not positive semi-definite. However, their slight 
modifications given by 

k{x, y) := ip{x -y)- ip{y) - ip{x) + ip{0) 

define positive semi-definite kernels, and the corresponding RKHSs were characterized 
in [16, Theorem 10.18 ]. By Remark [H Ex and Ej^ have the same minimizers, so that 
we can work with the original kernel K in the energy functional. 

2. Let Q, := T^. Up to an additive constant, Wahba's spline kernels are given by 
Kix, y) := '-^^B,^{\x " 2/1) = 2 ^ cos(2vrfc(x - y)) 

where i?2m. denotes the Bernoulli polynomial of degree 2m. Note that -B2m(l — i) = 
B2m{t). For example, we have 

B^{t) = t' -t + l, B^{t)=t^-2t^ + t^-^. 

Let us now introduce a slight generalization of ([5]), where we consider different kernels for the 
attractive interaction (attraction of the particles by the image profile w) and the repulsive 
interaction (particle-particle repulsion). In particular, we introduce the functions (p and ip 
related to different radial kernels and the generalized energy functional 




The following remark gives an intuition on the behavior of the corresponding quantization 
process. 

Remark 2 Assume that the kernel K in ^ is in addition continuous and an element of 
L^(ri X Q.) (Mercer kernel). Then it can be expanded into an absolutely and uniformly conver- 
gent series, 

oo 

K{x,y) = Kx{x,y) := ^eVe{x)r]i{y) 

of orthonormal eigenf unctions rji S L^(il) and associated eigenvalues Xi > of the compact, 
self-adjoint integral operator Tk on Lp'{yi) given by 

TKf{x) := / K{x,y)f{y) dy. 
Jn 
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Assume further that w can be also expanded into an absolutely convergent series w{x) 
Yl'kLi'f^k'nkix). Then the functional ([3]) becomes 



N „ . N 



Ek{p) = y2 [ w{x)KxiPi,x) dx - ^ S2 Kx{pi,Pj) 
i=i ^ ^,j=l 

N „ CO oo \ ^ 

= XI / ^^kVkix) ^ Xir]i{pi)r]i{x) dx - - ^ KxiPi,Pj) 

i=l k=l £=l i,j=l 

N oo . ^ N 

= X] X] '^k>^^^dPi) / r]k{x)rit{x) dx - - ^ Kx{pi,Pj) 

i=l k/=l ''^ i,j=l 

N oo ^ N 

= "^"^WkXkVkiPi) - 2 ^ Kx(.Pi,Pj) ■ 

i=l k=l «J=1 

On the other hand, if we consider Ek for another function {^{x) := Yl^=i'^kVk'nk{x) , we have 

N ^ > Af 



ek{p) = I '^ix)K\iPi^x) d^;- ^ KxiPi,Pj 

TV oo A ^ 

= ^Y^^^WkVkXk'nkiPi) ~ 2 ^ Kx{pi,Pj) 

i=l k=l i,j=^ 
N ^ N 

= y2 w{x)K^{pi,x) dx- - V Kx{pi,Pj) 
i=i 2 ^ .^^ 

with Kfj^{x,y) := X^^^i and fik '■= Vk^k; where we assume absolute convergence of 
the involved series. Hence, using a smoother kernel for the interaction with the datum w 
than Kx (i.e., fik decays faster than Xk) leads to the approximation of a smoother function w 
(wkVk decays faster than Wk), and vice versa. 



3 Mean-Field Limit 

We are interested in the passage to the hmit when the number of particles tends to infinity. 
For simphcity, we restrict our attention to radial kernels and $7 = M'^ although the analysis 
works as well for the periodic settings O = T and = T'^, d > 2 with the tensor product of 
radial kernels. Moreover, without loss of generality, we prescribe the normalization 

/ w{x)dx = 1, NX = 1 
Jn 

and suppose that > is compactly supported. 
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3.1 Passage to the Mean-Field Limit 

The evolution of the A'^-particle system according to the gradient flow of the discrete generalized 
energy functional 

N ^ N 

E{p) = E^^^{p) := X] / -x)dx- — ^ V(Pfc - Pi) (7) 

fc=l'^^ k,l=l 

is given, under the assumption that ip,^^ € C^{il), by 

1 ^ 

u;(x)V(/j(p,(t) -x)dx + -^VV(p.(t) -p^t)), i = l,...,N, (8) 

e=i 

subject to the initial condition 

Pi(0)=p°, i = l,...,N. (9) 

The mean field limit is obtained as the number of particles tends to infinity. Then, the vector 
of time-dependent particle positions p{t) € is replaced by the time-dependent probability 
measure f{x, t), where, roughly speaking, f{x, t) dx can be understood as the probability that 
a particle is located in the space element dx around the position x € at time t > 0. 

In the following, let A4{Q) denote the space of Radon measures on and Cd^) the space of 
continuous, compactly supported functions on fi. Further, let L°°(M_|_, denote the space 

of functions from to which are essentially bounded, i.e., f : t ^ /('lO — ft with 

esssupjeR^ d\f{-,t)\ < oo. Note that L°°(R+,M{n)) is the dual space of L^(R+,Cc{^)), 
the space of functions from to Cc{^) such that \\gi-,t)\\oo dt < oo, see, e.g., [2]. More- 
over, let us denote by the set of probability measures on Q,, i.e., / G A4^{Q) if and 
only if / is a nonnegative Radon measure such that df = 1. 

For any G N, let us denote by the empirical measures 

1 ^ 

1=1 

corresponding to the evolution of the A^-particle system ([5]). Then, each f'^ is a time-dependent 
probability measure, such that € L°°(]R+, A^(r2)). In the following theorem we carry out 
the rigorous mean field limit passage N ^ oo. 

Theorem 1 Let ipjip ^ C^{Q), where SJip is in addition bounded. Let {f^)N&N be given 
by (jlOp . corresponding to the system ([8]) with the initial datum ([9]). Moreover, assume that 
there exists a probability measure fo G such that /^(•,0) — > /o(-) weakly-* in as 

N ^ oo. 

Then there exists a subsequence which converges weakly-* in L°°{W^, 

to a time- dependent probability measure f G L°°{^^,Ai^{0,)) which solves, in the sense of 
distributions, the mean-field equation 

dtf = Vy( f {w{x)V^{y-x)- f{x,t)V^{y-x))f{y,t)dx 

= V-(V/C[/]/), (11) 
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where 

subject to the initial condition 



{w{x)ip{y -x)- f{x, t)^{y - x)) dx , 



/(•,0) = /o. 



(12) 



(13) 



Proof: First we show that for ah € N the empirical measures (jlOp are distributional 
solutions of ([IT]) . Indeed, considering a smooth, compactly supported test function ^ € C^(Ox 
[0, oo)), we obtain 



Jn 



1 

iV 



/^(y,i)5t^(y,i)dydt+ / /^(y,0)^(y,0)dy 

Jn 

N ^ N 

Y^^ti{p,{t),t)dt + -Y,^{p^{^),^) 



i=l 

N 



i=l 



i=l 

N 



N 



1 1''^ d d 1 1 " 

nL ^ ^e(PiW,i)-Vp,e(Pi(i),i)-^Pi(i) dt + -Y,^{MO),0) 



wix)Vip{piit)-x)dx--Y, 

Jo j^]^ yn 

/ / /^(y,t)V,^(y,t)- 
Jo Jn 



i=l 
N 

' ViP{pi{t) - pi{t)) 



dt 



w{x)V(p{y — x) dx — / / {x,t)Vip{y — x) dx 



Therefore, we have the identity 

f''{y,t)dtay,t)dydt+ [ /^(y,0)C(y,0)dy 

Jn 



Jn 

oo 



Jn 



f''iy,t)Vyay,t) 



w{x)V ip{y — x) dx — I {x,t)Vip{y — x) dx 



dy dt, 



(14) 



dy dt, 



for all test functions ^ G C^{Q x [0, oo)), that is the distributional formulation of (jlip with 
in place of / and the initial condition /^(•,0) = Xli^i ^i' ~ Pi)- 

Now, since {f'^)NGN is a sequence of time-dependent probability measures, it is uniformly 
bounded in L°°(R+, A^(r2)), so that there exists a subsequence {f^'')k£N which converges 
weakly-* to some / in L°°(M+, 7M^(r2)). We show that / is a distributional solution of ([TT]) . 
The limit passage in the linear terms of (jl4p follows immediately. Moreover, since Vip is as- 
sumed to be continuous and bounded, the sequence {x, t)V'il;{y — x) dx is uniformly equi- 
continuous and uniformly bounded for y € supp ^{■,t) and for almost every t £ M+. Therefore, 
due to the Arzela-Ascoli theorem, this sequence converges strongly in L°°(supp i{-,t)) for al- 
most all t G M_|_. Finally, the bounded convergence theorem ensures the strong convergence of 
the sequence in L^(R_|_, L°°(supp and this justifies the limit passage in the nonlinear term. 

Therefore, in the limit ^ oo, we have obtained 



Jn 

oo 



f{y,t)dti{y,t)dydt+ / /o(y)e(y, 0) dy 



(15) 



/o Jn 



f{y,t)Vy^{y,t) 



w{x)Vip{y — x) dx — / f {x,t)Vil}{y — x) dx 



dy dt . 
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which is the distributional formulation of (jlip subject to the initial condition (jl3p . ■ 

Equation pT]) describes the evolution of the time-dependent probability measure f{-,t) 
due to the mutual repulsive interaction between the particles and the attractive interaction 
with the datum w. In the following lemma we show that, under mild regularity assumptions 
on if, ip and /o, the solution / of pT]) is in fact classical. 



Lemma 1 Let V9J and S/tp be globally Lipschitz continuous on W^, i.e., there exist constants 
Li, L2 such that 

\Vip{x) -Vip{y)\ < Li\x-y\, 
iVV'(x)- VV'(y)| < L2\x-y\, 

for all X, y G ri. Let fo G C^{Q,) be nonnegative, compactly supported and fulfill fo{x) dx = 1. 
Then the corresponding distributional solution f of (jlip subject to the initial condition ()13p 
is in fact a classical solution with f e C^{Q x M+) and f{-,t) > for all t > 0. Moreover, 
f{-,t) is compactly supported on Q for any t € M+. 

Proof: Since the distributional solution / G L°°(M+, A^^(il)) constructed in Theorem [1] is a 
time-dependent probability measure, we have f{-,t) > and 

/ df{t,-) = l forant>0. 

The essential point is to observe that due to the assumptions on ip and ip, the transport field 
V/C[/] for (fTT]l is Lipschitz continuous. Indeed, we have 

\VIC[f]{p) -VJC[f]{q)\ < [ w{x)\Vip{p-x) -V^{q-x)\dx 

+ / f{x,t)\\/'ip{p — x) — \/ilj{q — x)\dx 
Jn 

< Li w{x)\p — q\dx + L2 / f{x,t)\p — q\dx 

Jn Jn 

< {Li + L2)\p - q\ . 

Therefore, / is a solution of a hyperbolic transport equation with globally Lipschitz-continuous 
transport field V/C[/]. As such, the values of the initial condition /o propagate along the 
characteristics 

x(t) = V/C[/](x(t)) (16) 

with finite speeds. Due to the assumption /q G C^{0,), from the standard theory of hyperbolic 
transport equations (method of characteristics and Cauchy-Lipschitz theorem, see [3]) it follows 
that / S C^(il X The compactness of the support of / for all times follows from the 

assumed compactness of the support of fo and the finite characteristic speeds (fT6]) . ■ 



Remark 3 Similarly as in Remark\^ we observe that the evolution induced by (jlip with two 
different interaction kernels </? and ip is in fact equivalent to an evolution produced by using 
same interaction kernels, but with a modified data w. Indeed, assuming that (p and ip are 
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continuously differentiable and that the Fourier transform of tp is nonzero almost everywhere 
(as is the case for positive definite kernels), the Fourier transform of (jlip reads 



dtf = i^ ■ {[uN^ - fVi^) * f) 

= -icp(v^(-| -/)*/). 

Applying the inverse Fourier transform, we get 

dtf = V- {{w - f) * W)/) , 

where w is the inverse Fourier transform of {w(p/ip). Therefore, taking ip smoother than ip 
corresponds to a smoothing of w, while ^ smoother than ip corresponds to a "sharpening" 
(anti- smoothing) ofw. 



3.2 Energy dissipation, long time behavior and equilibria 

Let us observe that the formal hmit of the discrete energy ([7]) as — > oo is given by the 
continuous energy functional 

£[f]= I I w{x)(f{p - x)f{p)dxdp-}- I I f{x)ip{p- x)f{p)dpdx , (17) 
Jn Jn ^ Jn Jn 

defined for all / G I^iVC). Let us mention that the corresponding formal gradient flow with 
respect to the topology of 2-Wasserstein distance on the space of probability measures (see ^lOj 
or [H] for details) is given by the hyperbolic transport equation (fTT]l . This suggests that the 
energy £'[/(x, •)] actually is a Lyapunov functional, thus nonincreasing along the solutions 
of (HIl): 

^^[/(x,t)]<0. (18) 
Indeed, at least for classical solutions, this inequality can be proven rigorously: 

Lemma 2 Let -0 &e of the form Let f G C^(r2x]R_|_) he a classical solution of pT]) - (fT3]) 
in the sense of Lemma [21 Then (llSp holds true. 
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4-Afix,t)] = I I w{xMp-x)Vp-{V}C[f]{p,t)f{p,t)) dpdx 



Proof: The proof follows from the direct calculation 
d 

Jn Jn 
1 

~2 
1 

2 Jn Jn 



ij{p - x)f{x, t)Vp ■ (V/C[/] {p, t)f{p, t)) dp dx 

n Jn 

4^{p-x)f{p,t)V^ ■ {V)C[f]{x,t)f{x,t)) dpdx 
w{x)Vip{p - x) ■ {WIC[f]{p,t)f{p,t)) dpdx 



n Jn 

+\ I I V^{p-x)f{x,t)-{VIC[f]{p,t)f{p,t))dpdx 
^ Jn Jn 

-\ \ \ VV(p-x)/(p,t)-(V/C[/](x,t)/(x,t)) dpdx 
^ Jn Jn 

w{x)Vip{p - x) -V^{p- x)f{x,t)) ■ {VIC[f]ip,t)fip,t)) dpdx 

VlC[f]{p,t)f{p,t)\^f{p,t)dp<0, (19) 

where we have used the symmetry of ip, i.e., Vip{p — x) = —Vtp{x — p), and integration by 
parts, where the boundary terms vanish due to the compact support of /. (In the case Q = T'^ 
due to the periodicity of the interaction kernels.) ■ 

In the rest of this section, we study the question whether the classical solution / of (lllh 
tends for t — )• oo to an equilibrium /*, characterized by the condition |V/C[/*](p)/*(p)|^ f*ip) = 
a.e. on fi, which stems from setting the right-hand side of (jl9p equal to zero. Without loss 
of generality (see Remark [3]), we restrict our attention to the case ip = ip] then we have 

V/C[/](p,t)= f {w{x)- f{x,t))V^{p-x)dx = {{w- f{;t))*V^){p). (20) 
Jn 

Moreover, we adopt the assumption that / is a classical solution of (llip in the sense of 
Lemma [TJ 

First we show that the energy functional (117p is bounded from below. 



Lemma 3 Let $ : [0, oo) — t- M 6e concave and monotone increasing and '5(0) he finite. Then 
E[f] with Lp = ip = ■ \) is bounded from below for all f G JvO'{^). 

Proof: Without loss of generality, we can assume that <^(0) = 0. Concavity and monotonicity 
of $ imply its subadditivity [H], which implies further 



^{\p - x\) < <^{\p\) + $(|x|) , $(|p - x|) > $(|p|) - ^-dx] 
Therefore, since w, f >0 and w{p) dp = Jq f{p)dp = 1, we obtain 



£[f] > 



n Jn 



w{xm\p\) - $(1x1)) - \f{xm\p\) + <i>(ix|)) 



f{p) dpdx 



7i;(x)<I>(|x|) dx , 
which is the announced boundedness from below. 
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Remark 4 Note that if^{r) grows as r'^ with r > 2 as r — >■ cxd and w is compactly supported, 
then £ is not hounded from below. This can be seen as follows: Define fq{x) := ^{d{x — q) + 
S{x + q)) for q G M*^. Then, due to the compact support of w, we get 

Sif,] = ljwix)mq-x\+^\q + x\2))dx-^m<l\) 

~ i^r - \{Mr = (1 - ^2oigr ^ -oo as ^ oo . 

Now, since by Lemma [2] the energy £[f{x, t)] is nonincreasing as a function of time and by 
Lemma El bounded below, the hmit of £[f{x,t)] as t — ?> oo exists and is finite. 

Due to the boundedness of f{-,t) in the space of Radon measures, there exists a sequence 
tj — > oo and a Radon measure /* G such that f{-,tj) — )• /* weakly-* as tj oo. Since, 

as assumed, / is a classical solution with / G C^{Q, x M+), we also have £[f{-,tj)] — )• £[f*] 
as tj — >■ oo, and, consequently, £[f{-,t)] — > £[f*] as tj — > oo. By ([T9|), the equilibrium /* is 
characterized by the condition 

|V/C[r](p)r(p)|V*(p) = a.e. onQ. (21) 

By ()20p the choice f* = wis always a solution of ()2ip . This corresponds to the intuitive 
expectation that, as we let the number of particles tend to infinity, we should recover the 
profile w in the long-time limit, regardless of the initial distribution of particles. However, the 
question whether the choice f* = w is the unique solution of (f2T]l in the class of probability 
measures seems to be rather nontrivial. Although we believe that the affirmative indeed holds 
for a broad class of interaction potentials, we are so far only able to provide a proof for 
the special case = M and (p{-) = ip{-) = | • |. Unfortunately, this case does not match 
our assumptions on the smoothness of (p and ^ made in Theorem [T] and Lemmata [T] and [21 
However, the weak formulation (jlSp of equation (jlip perfectly makes sense if we insert the 
distributional derivative <p'{-) = ip'{-) = sign(-) into (fT2]l . as long as w,/ G L^{Q,), since then 
the integrals 

/ w{x)ip' {p — x) dx = / w{x)sign{p — x) dx 
Jn Jn 

and 

f {x)iIj' {p — x) dx = / /(x)sign(p — x) dx 
n Jn 

are well defined for all p G M and uniformly bounded. Consequently, we can formulate the 
following Lemma: 

Lemma 4 Let = M and (p{-) = = | • [. Let w > be compactly supported in Q 
and such that J^w{x) dx = 1. Then the solution f* = w of (|2ip is unique in the class 
X ■.= {f e L^{n) : / > a.e. on n,J^f{x)dx = l} . 



Proof: Let f £ X fulfill (|2ip . which can be recast as 



n 



{w{x) — /(x))sign(p — x) dx 



f{p) = a.e. on (22) 
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Due to the normalization w{x) dx = f{x) dx = 1, we have 




/ /(x)dx+ / /(x)dx 



= 2W{p) - 2F{p) , 



where we denoted 



■p 



W{p) : 



w^x) dx 



and F{p) := / /(x)dx 



— oo 



— oo 



The condition ()22p can then be rewritten as G{p)f{p) = for almost all p G 0; let us note 
that G is a continuous function. By assumption, there exist real numbers a and uj such that 
supp w C [a,a;]. We prove that G = in three steps: 

• First we show that G{p) = for all p < a. For a contradiction, let us assume that 
there exists a po < a such that G{po) 7^ 0. Then we have G{po) = 2W{po) — 2F{pq) = 
—2F{pq) < and, therefore, there exists a set S C {—oo,po) of positive Lebesgue 
measure, such that / > almost everywhere on S. But then F > almost everywhere 
on S, and, consequently, G = —2F < almost everywhere on S, a contradiction to 



for all p > UJ. 

• Finally, we prove that G{p) = also for all p G {a,uj). By the continuity of G, the 
set S := {p G (a,w) : G{p) ^ 0} is open. Since (j22]) dictates that /(p) = almost 
everywhere on 5, we have for every po S and 6 > small enough. 



Therefore, G is nondecreasing on S, and, thus, G is nondecreasing everywhere on (q,w). 
Since G{a) = G{uj) = and G is continuous, we conclude that G = 0. 

We finish the proof by observing that 



G{p)f{p) = 0. 



• Using a symmetry argument, we can prove that 






for all p G M 



implies f = w almost everywhere on M. ■ 
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4 Numerical Examples 



In this section we present few numerical results for the discrete particle system ([6]) and the 
mean- field limit (jlip . We consider two cases: 

• Smoothing case: ip{s) = \s\^'^ and Tp{s) = \s\ 

• Sharpening case: ip{s) = \s\ and ip{s) = \s\^'^ 

In both cases, we consider the ID full-space setting = M with the datum w = 4x[q 25,5]- 
For the discrete particle system, we use = 20, 50 and 100 particles. We integrate the ODE 
system ^ in time using the explicit Euler method until the steady state (which does not 
depend on the initial condition). The results for the smoothing and sharpening cases are 
shown in Fig. [2j 



0.1 0,2 3 0.4 0.5 0.6 0,7 0.25 0.3 35 4 0.45 0.5 




Figure 2: Steady states of the discrete particle system ([8]) with N = 20, 50, 100 particles (first, 
second and, resp., third row) in the smoothing (left panels) and sharpening case (right panels). 
The datum w is visualized with the dashed line. Note the different horizontal axis limits in 
the smoothing and sharpening case. 

For the mean-field limit ([TT]) . we impose the initial condition /o = 4x[o. 65,0.9]- We dis- 
cretize (jlip using the semi-implicit finite difference method with upwinding in space and ex- 
plicit Euler method in time. Snapshots of the solutions are shown in Fig. [3] for the smoothing 
case and in Fig. [5] for the sharpening case. 

As we expect, the solutions converge to some steady states as t — > cx) in both the discrete 
and mean-field cases. In the smoothing case, the equilibrium profile is a smoothed version of 
the data w, while in the sharpening case the equilibrium is an anti-smoothed version of w. It 
is interesting to compare these results with the discrete particle calculation from the previous 
example. Indeed, the steady state particle distributions shown on the left panels of Fig. [2] can 
be regarded as approximations of the steady state on Fig. [3l and the same hold for the right 
panels of Fig. [2] and the steady state of Fig. [H 
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0.2 0.4 0.6 0.8 1 0.2 0.4 0.6 0.8 1 



t = 0.0 t = 1.5 




0.2 0.4 0.6 0.8 1 0.2 0.4 0.6 0.8 1 



t = 3.0 t = 5.0 

Figure 3: Smoothing case of the mean-field Hmit The sohd fine represents the solution 

/, the dashed line the data w. The upper left panel shows the initial condition, the lower right 
panel is the steady state. 

In our last example, we show how the kinetic equation (jlip can be used for deblurring 
of images in the case when the blurring kernel is known. As our "image" we take again the 
datum w = 4x[o.25,5]) and blur it by applying 200 time-steps of ([TT]) with ip{s) = \s\^'^ and 
ijjis) = \s\ (smoothing case), the time-step length is 0.01. For simplicity, we use the initial 
datum /o := w. The result of this blurring process, g = f{-,t = 2), is shown on Fig. [5l left 
panel. Then, we perform deblurring of the "image" by applying (|lip with reversed interaction 
potentials, i.e., 99(5) = |s| and ip{s) = \s\^'^, and with w := g. Again, for simplicity, we take 
fo := w = g as the initial condition for /. We let (jlip evolve until steady state, which is 
reached around t = 30, and plot it in the right panel of Fig. [5l We see that f{-,t = 30) is 
indistinguishable from the original, un-blurred image w; the relative difference between them 
in the L^-norm is approximately 0.6%. 

Let us mention that numerical implementation of (jlip in the spatially 2D setting, which 
potentially might be of interest for application in image processing (deblurring of images), is 
quite a demanding task. The main reason is the high numerical cost, caused by the necessity 
of evaluation of the convolution / * V'ip in each time step, which in general takes 0{N'^) 
multiplications if the grid consists of N'^ points. A possible speed-up of this operation can be 
achieved by using fast multipole expansion or FFT-based methods, see for instance [5]. We 
postpone this task for future work. 
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Figure 4: Sharpening case of the mean- field limit The solid line represents the solution 

/, the dashed line the data w. The upper left panel shows the initial condition, the lower right 
panel is the steady state. 
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